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Abstract. 

f I I ' We examine structural relaxation in a supercooled glass-forming liquid 

O simulated by NVE molecular dynamics. Time correlations of the total kinetic 

ty^ , energy fluctuations are used as a comprehensive measure of the system's approach 

to the ergodic equilibrium. We find that, under cooling, the total structural 
C^ ' relaxation becomes delayed as compared with the decay of the component of the 

e' intermediate scattering function corresponding to the main peak of the structure 

factor. This observation can be explained by collective movements of particles 
T-i—i , preserving many-body structural correlations within compact 3D cooperatively 

_H 1 rearranging regions. 
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Understanding the microscopic mechanisms of the super- Arrhenius slowing down 
in fragile glass-forming liquids [T] remains a major problem in condensed-matter 
physics [5] [3]. It is understood that below the temperature T4 that marks the crossover 
to super-Arrhenius behaviour liquid dynamics is controlled by the topography of 
the potential energy-landscape (PEL) [4]. The PEL is known to be divided into 
metabasins - areas of PEL confining sets of mutually well-connected local energy 
minima [5]. Inter- metabasin transitions, mediating the (primary) a— relaxation, 
involve uncorrelated movements of large groups of particles [6]. On a shorter time 
scale, confinement to a metabasin imposes structural constraints upon the atomic 
motions, reducing the number of accessible degrees of freedom. This renders intra- 
metabasin dynamics highly collective, in accordance with the Adam-Gibbs concept of 
independent cooperatively rearranging regions (CRR) 0. 

The geometry, and dynamics of CRR remain elusive. By definition [7], a CRR 
is a minimum-size region of structure possessing a relaxational degree of freedom 
independent of the region's environment. It was suggested [S] that this degree of 
freedom can be associated with a coherent linear translation of a string-like cluster 
of particles within slowly-changing structural environment [9] . On the other hand, 
the random first-order transition theory of glasses [10] predicts a transition, at some 
stage of the liquid's cooling, from a string-like shape of CRR to a compact one. 
The concept of a compact CRR presumes that there exists a relaxational degree of 
freedom that makes it possible for a 3D configuration of particles confined to that 
CRR to perform a cooperative movement that would preserve some of the many-body 
structural correlations constraining the configuration. 

A distinctive aspect of the relaxation dynamics associated with compact CRR is 
that the described cooperative 3D movements of particles are expected to result in 
decorrelation of the density fluctuations while preserving slower-decaying many-body 
structural correlations. The relaxation of the density fluctuations is described by the 
intermediate scattering function F{Q,t)^^. Its slowest-decaying component (within 
the relevant range of Q) corresponds to Qm - the position of the main maximum of 
the structure factor S{Q). In the case of cooperative 3D particle movements within 
compact CRR, F{Qjn,t) is expected to decay well before the total liquid's structural 
relaxation has been accomplished. This can be understood using the following simple 
model. Consider slowly relaxing structural domains (clusters) immersed in a faster 
relaxing liquid. Let the clusters perform uncorrelated rotational and translational 
movements. It is easy to see that F{Qm,t) for such a model would decay well before 
the clusters dissipate as structural entities. 

Detecting the described effect requires the use of a comprehensive measure of 
the structural relaxation that would be independent of the dissipation of density 
fluctuations. In this paper, we report a molecular dynamics (MD) simulation of a 
fragile glass-forming liquid where the structural relaxation is measured by the time 
correlations of the fluctuations of the system's total kinetic energy. We find that, as 
the hquid is cooled sufSciently close to the mode-couphng theory (MCT) [12] critical 
temperature Tc, its total structural relaxation becomes delayed as compared with the 
respective decay of F{Qm, t). This anomaly can be interpreted as an indication of the 
existence of slowly dissipating many-body correlations constraining the movements of 
cooperatively rearranging 3D groups of particles. These groups can be identified with 
compact CRR. 

The MD simulation we report here explores a fragile simple one-component glass- 
forming liquid demonstrating a pronounced tendency for icosahedral clustering [13] 
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Figure 1. Asymptotic behaviour of Q(t) as predicted by|4] Solid line, T = 0.75; 
dashed line, T = 1; chain-dotted lino, T = 1.35. Inset: solid line, r as evaluated 
from the asymptotic slope of Q(0)/Q{t), defined in |4] boxes, r calculated by 
integrating c{t), defined in[5l 



(named Z2 in that reference), with the estimated T^. ~ 0.65, fragihty index B — 4.5, 
and T4 = 1.1. For each temperature, an equihbrium simulation in an NVE ensemble 
was performed using a system of 16384 particles. A system of 128000 particles was 
also tried, and no size-dependent effects have been observed. 

Consider a system of N particles in an NVE ensemble. Using the expression for 
the isometric heat capacity |14l 115) , its entropy s can be linked to the variance of the 
kinetic energy k (both per particle): 

-1 - M -M '^ 2 (P) - (fc)^ m 
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The respective expression for the entropy explored by the system within the time 
interval t, s{t), [TH], will be: 
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where ( )t denotes averaging over the time interval t. The system's approach to the 
ergodic equilibrium can be comprehensively assessed as the difference between the two 
quantities [17]: 
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Note that Q.{t) is a formal analogue of an earlier suggested measure of ergodic 
convergence [18] calculated from the time-dependent variance of the single particle 
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Figure 2. Evolution of the correlation function c(f), as defined by [6] shown 
as a function of the mean-square displacement. Solid line, T = 1; dashed line, 
T = 0.75; chain-dotted line, T = 1.35. Inset: the short-time behaviour of the 
same curves. Note that in the inset they are shown as a function of time. 



energy. Following the arguments presented in that work, it is also possible to conclude 
that the large-i asymptotic behaviour of ri(t) can be described as 

(k)' ml) 



n{t)/n{Q) 



(fc}2 - (fc2) 
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where 



dic(i)/c(0) 
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and c{t) is the time correlation function for the system's total kinetic energy 
fluctuations: 

c(i) = {5k{t)5k{Q)) (6) 

with 5k{t) — k{t) — (k). This asymptotic behaviour of ri(t) is indeed observed in 
our MD simulations, both above and below Ta, see Fig. 1. The inset of that figure 
demonstrates the agreement between r as evaluated from the asymptotic slope of 
ri(0)/ri(i) m and that obtained by direct integration of the correlation function c{t) 

m 

Thus, the system's relaxation, as measured by fi{t), is ultimately controlled by the 
decay of c{t). The general structure of this correlation function is shown in Fig. 2. Its 
short-time behaviour (shown in the inset) is dominated by large-scale oscillations. The 
latter are evidently caused by particles vibrations within the cages of their immediate 
neighbours, the frequency being determined by the mean-square force J19J with weak 
temperature dependence. 

The described short-time behaviour of c{t) obviously excludes the possibility of 
using its normalised integral r, equation (5), which determines the decay rate of 
fl{t), as an unambiguous measure of the liquid's relaxation rate. However, the initial 
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Figure 3. Structural relaxation in the temperature domain of supercooled 
liquid dynamics. Solid lines, c(t); dashed lines, F{Qm,t); left and right curves, 
respectively, correspond to T = 0.75 and T = 0.7. Dots: the fits of KWW 
stretched exponential approximation Aexp [— (i/r)'']. The values of the fitting 
parameters are presented in Table 1. All sets of data presented in the plot are 
scaled by the respective values of A. 



oscillatory regime of c{t) is distinct in time, and decoupled from the relaxational regime 
that follows it. The structural relaxation in the asymptotic a-regime can therefore 
be adequately measured by the long-time decay of c(i). We remark that, above Ta, 
this longer-time relaxational evolution of c(i), when presented as a function of the 
mean-square displacement, Fig. 2, demonstrates the standard liquid universality of 
the diffusion-relaxation relation. This universality evidently breaks below Ta- 

We now exploit c(t) for a detailed analysis of the a-relaxation process in the 
supercooled dynamics regime. In Fig. 3, the asymptotic behaviour of c(t) is compared 
with that of F{Qm,t), for two temperatures below Ta- For convenience, both are 
plotted as a function of the mean-square displacement. To facilitate the comparative 
quantitative analysis, the correlation functions presented in Fig. 3 have been fitted, 
within the relevant domains of time, with the Kohlrausch- Williams- Watts (KWW) 
stretched exponential function Aexp[—{t/T)^] which we use here in the equivalent 
form Aexp{ — [{r-^{t))/6DT]^}, D being the diffusion coefficient. The parameters' 
values for the fits are compiled in Table 1. Furthermore, to make it possible to 
compare the two kinds of correlation functions on equal footing, each function has 
been scaled by the respective value of the KWW parameter A. Both F{Qm,t) and 
c{t) appear to be well described by the KWW approximation within the a-relaxation 
time domain. The latter follows the initial stage of the relaxation process, confined 
to {r'^{t)) < 1, that corresponds to conformational rearrangements of particles within 
the cage of the first neighbours. 

For T = 0.75, the asymptotic behaviour of the two correlation functions is 
indistinguishable, within the accuracy of the KWW fit [20]. However, as the liquid 
is cooled to T = 0.7, the asymptotic decay of c{t) becomes considerably delayed as 
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Table 1. Parameters' values for the stretched-exponcntial KWW fits to the 
relaxation functions shown in Fig. 3 
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Figure 4. Top panel: radial distribution function for the simulated liquid at 
T = 0.7. Bottom panel: F(Q,t). From left to right: Q = 4, Q = 6, Q = 2, 
Q = 0.6, Q = 7.2. Dots correspond to the fits of KWW. Inset shows the KWW 
parameters derived from the fitting. Open circles: r (left-hand scale); dots: /3 
(right-hand scale). 



compared with the respective decay of F{Qm,t). We also note that the observed 
delayed decay of c(t) is accounted for by a measurable reduction of the respective 
KWW stretching parameter /3 (Table 1). The relaxation stretching, arising under 
cooling below Ta, is commonly viewed as a result of superposition of relaxation 
processes with an extended range of relaxation times. Thus, the additional relaxation 
stretching detected using c(t) indicates an additional slow relaxation process not 
captured by F{Qm,t). 

We remind that a system with a finite configurational entropy possesses a 
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finite range of structural correlations. Therefore, the ergodicity-restoring structural 
relaxation in a liquid is a strictly local process. In the case of pair correlations, 
the correlation range r^ corresponds to vanishing radial distribution function g{r). 
Respectively, if relaxation of the pair correlations is considered in terms of F{Q,t), 
the lowest Q relevant for the structural relaxation can be estimated as Qmin = "^n/rc- 
Q-dependent relaxation time of F{Q,t) can be estimated from the structure factor 
S{Q)^ as t{Q) = S{Q)/Q'^. In a typical dense simple liquid llj S{Q,n)/S{Qmin) > 
{Qm/Qmin)^, and, therefore, T{Qm) is expected to exceed T{Qmm)- In Fig. 4, we 
test this conjecture for the simulated liquid at T = 0.7. The figure shows that g{r) 
apparently vanishes beyond re = 10, from which we estimate Qmin ~ 0.6. Relaxation 
of F{Q,t) is analysed within the relevant range of Q bounded by Qmin and Qm- For 
each value of Q, the decay of F{Q,t) is quantified by the KWW fit; the Q variations 
of the resulting KWW parameters are shown in the inset. Evidently, T{Qm) is a 
pronounced absolute maximum within the explored range of Q, whereas f3 doesn't 
change significantly. Therefore, structural relaxation in our liquid at T = 0.7, in the 
pair approximation, is entirely controlled by F(Qm,t). This leads to the conclusion 
that the observed delay of the total structural relaxation relative to the decay of 
F{Qm,t) is due to the delay of the relaxation of higher-order correlations relative to 
that of pair correlations. 

The main peak of S{Q) is linked to the local order (cages of nearest neighbours). 
Therefore, one reason for the decay of F{Qm,t) can be a comprehensive dissipation 
of the local order due to uncorrelated particle motions. But it can also be caused 
by collective particle movements preserving the higher-order structural correlations. 
This will delay the total structural relaxation relative to the decay of F{Qm, t) as we 
observe for T = 0.7. The increase in the life-time of these structural correlations can 
be attributed to a rapid buildup of locally preferred structure that is expected upon 
cooling the liquid towards Tc- This transformation in the relaxation dynamics can be 
interpreted as a crossover to the regime of compact CRR as conjectured in Ref. [TO] . 

The locally preferred structure in this system is icosahedral. Like an earlier 
studied liquid with a similar structure [551 123]) it demonstrates a strong tendency 
for a low-dimensional icosahedral aggregation [T3] growing towards percolation as 
T — » Tc- Relative mobility of these structural elements can be conjectured as a 
conceivable reason for the decoupling between c(t) and F{Qm,t). As a corroborating 
observation, the anomalous decay of F{Qm,t) reported here can be compared with a 
related Q-dependent anomaly in the non-ergodicity parameter in a supercooled low- 
density liquid approaching gelation transition [53] caused by a percolating cluster 
network. Compact 3D CRR were also found below T^ in the binary Lennard- Jones 
(BLJ) system [25] and NiZrj26]. and bond-preserving movements of structural domains 
was observed at a saddle point crossing of the BLJ system [27] . On the other hand, a 
study of the frequency-dependent specific heat in silica [28] , also based on the analysis 
of the kinetic-energy fluctuations, showed no discrepancy between the relaxation of 
these fluctuations and the decay of the self part of F{Q, t). 

We note that, because of the particles' indistinguishability, the structural 
correlations can survive particles exchange induced by vacancy hopping. Therefore, 
the slow energy fluctuations we report here cannot be necessarily associated with 
breaking bonds connecting the nearest neighbours. The vacancy-assisted particle 
hopping, common in supercooled liquids, significantly complicates the problem of 
identifying the groups of particles involved in the conjectured compact CRR dynamics. 

One of the consequences of compact CRR concerns the relaxation of shear stress. 
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The liquid can possibly accommodate to the shear strain, and relax the induced 
stress by a mutual rearrangement of slowly relaxing clusters [H] before their internal 
structure has dissipated. Thus, shear stress relaxation in a supercooled liquid, like 
the dissipation of the local density fluctuations as measured by F{Qrmt), doesn't 
necessarily imply structural relaxation. 

In summary, we found a novel aspect of the supercooled liquid dynamics 
manifested by the advanced decay of F{Qm,t) relative to the actual structural 
relaxation. The latter has been assessed from the decay of the time correlations of the 
system's total kinetic energy fluctuations. The observed effect can be interpreted as 
being caused by the cooperative movements of particles constrained by slowly-decaying 
many-body structural correlations within compact CRR which appear as the liquid is 
cooled close to Tc- The observed transformation in the relaxation mechanism can be 
regarded as a result of the metabasin topography of the PEL of our model. 

The authors thank Prof. D. J. Wales for very useful discussions. The support 
from the Centre for Parallel Computers (PDC) is gratefully acknowledged. 

References 

[1] C. A. Angell, J. Non-Cryst. Solids 131-133, 13 (1991). 

[2] M. D. Ediger, C. A. Angell, S. R. Nagel, J. Phys. Chem. 100, 13200 (1996). 

[3] P.G. Debenedetti, F.H. Stillingcr, Nature 410, 259 (2001). 

[4] M. Goldstein, J. Chem. Phys. 51, 3728 (1969) 

[5] F. H. Stillinger, Science, 267, 1935 (1995). 

[6] G.A.Appignanesi, J. A. Rodriguez Fris, R.A. Montani and W. Kob, Phys. Rev. Lett. 96, 057801 
(2006) 

[7] G. Adam and J.H. Gibbs, J. Chem. Phys. 43, 139 (1965) 

[8] N. Giovambattista, S.V. Buldyrev, F.W. Starr, and H.E. Stanley, Phys. Rev. Lett. 90, 085506 
(2003) 

[9] C. Donati, J.F. Douglas, W. Kob, S.J. Plimpton, P.H. Poole, S.C. Glotzer, Phys. Rev. Lett. 
80, 2338 (1998) 
[10] J.D. Stevenson, J.Schmalian, P.G. Wolynes, Nature Physics, 21, 268 (2006) 
[11] J. P. Hansen and I. McDonald, Theory of Simple Liquids, (Academic Press, London, 1976) 
[12] W. Gotze and L. Sjogren, Rep. Prog. Phys., 55, 241 (1992) 
[13] J.P.K. Doye, D.J. Wales, F.H.M. Zetterling, and M. Dzugutov, J. Chem. Phys. 118, 2792 

(2003) 
[14] J.L. Lebowitz, J.K. Percus and L. Verlet, Phys. Rev. 153, 250 (1967) 
[15] J.M. Haile, Molecular Dynamics Simulation: elementary methods, Wiley & Sons (1992) 
[16] S. K. Ma, Statistical mechanics, (World Scientific, Singapore, 1990). 
[17] R. G. Palmer, Adv. Phys. 31, 669 (1982). 

[18] R. D. Mountain, D. Thirumalai, J. Phys. Chem. 93, 6975 (1989). 
[19] J. P. Boon and S. Yip, Molecular Hydrodynamics, McGraw-Hill, New York (1980) 
[20] The statistical accuracy of the correlation function c(i) as shown in Fig. 4 required the averaging 

over the total simulation run of 10 timesteps. 
[21] P. G. de Gennes, Physica, 25, 825 (1958) 

[22] M. Dzugutov, S.I. Simdyankin, and F.H.M. Zetterling, Phys. Rev. Lett., 89, 195701 (2002) 
[23] Y. Gebremichael, M. Vogel, M.N.J. Bergroth, F.W. Starr, and S.C. Glotzer, J. Phys. Chem, 

109, 15068 (2005) 
[24] E. Zaccarelli, S.V. Buldyrev, E. La Nave, A.J. Moreno, L Saika-Voivod, F. Sciortino, and P. 

Tartaglia, Phys. Rev. Lett., 94, 218301 (2005) 
[25] K. VoUmayr-Lee and A. Zippelius, Phys. Rev. E, 72, 041507 (2005) 
[26] H. Teichler, Phys. Rev. E, 71, 031505 (2005) 

[27] S.Trygubenko and D.J.Wales, Journ. Chem. Phys, 120, 2082 (2004) 
[28] P. Scheidler, W. Kob, A. Latz, J. Hofbach, K. Binder, Phys. Rev. B, 63, 104204 (2001) 
[29] F.H StiUinger, Journ. Chem. Phys, 89, 6461 (1988) 



